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I.  INTRODUCTION 


This  project  has  involved  several  investigators  over  the  duration  of  the  con¬ 
tract  and  has  touched  on  several  fields  of  research.  The  overriding  theme  of  the 
research  has  been  on  grid  quality  and  the  improvement  of  numerical  simulation 
through  control  of  error.  This  has  been  achieved  by  making  major  modifications 
to  the  EAGLE  grid  code  and  the  EAGLE  flow  code  that  is  referred  to  as  the 
MISSE  code.  For  the  most  part,  these  modifications  are  transparent  to  the  user. 

In  the  earlier  form  of  the  adaptive  EAGLE  system  [1,2],  the  coupling  of 
the  adaptive  grid  system  with  a  computational  fluid  dynamics  (CFD)  code  re¬ 
quired  the  encapsulation  of  both  the  entire  grid  code  and  the  CFD  flow  code 
into  separate  subroutines  and  the  construction  of  a  driver  to  call  each.  This  was 
inefficient  in  that  it  included  some  unnecessary  parts  of  the  grid  code  and  re¬ 
quired  significant  modification,  and  perhaps  restructuring,  of  the  CFD  code.  In 
particular,  the  flow  code  arrays  and/or  the  grid  code  arrays  had  to  be  modified 
to  be  compatible  in  structure. 

In  the  effort  described  here,  the  elliptic  grid  generation  procedure  in  the 
EAGLE  grid  code  has  been  separated  from  the  main  code  into  a  subroutine, 
and  a  new  subroutine  that  evaluates  several  grid  quality  measures  at  each  grid 
point  has  been  added.  The  elliptic  grid  routine  can  now  be  called  either  by  a 
CFD  code  to  generate  a  new  adaptive  grid  based  on  flow  variables  and  quality 
measures  through  multiple  adaptation  or  by  the  EAGLE  main  code  to  generate 
a  grid  based  on  quality  measure  variables  through  static  adaptation.  Arrays  of 
flow  variables  can  be  read  into  the  EAGLE  grid  code  for  use  in  static  adaptation 
as  well. 
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These  major  changes  in  the  EAGLE  adaptive  grid  system  make  it  easier  to 

convert  any  CFD  code  that  operates  on  a  block-structured  grid  (or  single- block 

grid)  into  a  multiple  adaptive  code.  The  conversion  procedure  is  accomplished  by 

adding  the  elliptic  grid  generation  subroutine,  and  certain  other  subroutines  from 

the  EAGLE  grid  system  that  are  involved  in  the  elliptic  grid  generation  process, 

to  the  flow  code.  The  CFD  code  may  then  call  the  elliptic  grid  generation  routine 

at  each  time  step  when  a  new  grid  is  desired.  The  CFD  code  passes  its  current 

solution  to  this  EAGLE  routine  via  a  scratch  file,  and  the  EAGLE  routine  returns 

the  new  grid  to  the  CFD  code  also  via  the  scratch  file.  One  restriction  is  that  the 

initial  grid  must  be  generated  by  the  EAGLE  system  or  be  processed  through 

» 

that  system,  which  provides  the  necessary  parameters  and  structural  information 
to  be  read  from  files  by  the  adaptive  EAGLE  routines.  This  structure  eliminates 
the  need  for  compatibility  between  CFD  and  grid  arrays. 

In  the  present  work,  the  control  function  approach  is  used  as  the  basic 
mechanism  for  the  adaptive  grid  generation.  The  static  and  multiple  adaptive 
grid  generation  techniques  are  investigated  by  formulating  the  control  functions 
in  terms  of  either  grid  quality  measures,  the  flow  solution,  or  both. 

Previous  work  [1,2]  allowed  the  grid  to  only  adapt  to  the  gradient  of  a  vari¬ 
able.  The  work  described  here  has  extended  this  adaptive  mechanism  to  also 
allow  adaptation  to  the  curvature  of  a  variable  or  to  the  variable  itself.  The 
system  provides  for  different  weight  functions  in  each  coordinate  direction.  In 
addition,  the  mechanism  now  includes  the  ability  to  calculate  the  weight  func¬ 
tions  as  weighted  averages  of  weight  functions  from  several  flow  variables,  and/or 
quality  measures.  This  allows  the  adaptation  to  take  into  account  the  effect  of 
many  of  the  flow  variables  instead  of  just  one.  The  construction  of  the  weighted 
average  of  flow  variables  and  quality  measures  and  the  choice  of  adaptation  to 
gradient,  curvature,  or  variable  are  all  controlled  in  each  coordinate  direction 
through  input  parameters. 
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The  quality  measures  now  available  in  the  EAGLE  grid  system  are  skewness, 
aspect  ratio,  arc  length,  and  smoothness  of  the  grid.  These  grid  quality  measures, 
and  the  resulting  control  and  weight  function  values,  can  be  output  for  graphical 
contouring. 

In  addition  to  the  above  adaptation  using  the  elliptic  grid  equations,  a  grid 
refinement  capability  is  introduced.  When  a  grid  is  constructed,  the  MISSE  flow 
solver  may,  in  each  block,  compute  a  solution  on  the  input  grid  or  compute  a 
solution  on  a  refined  or  coarsened  grid.  This  grid  refinement  or  coarsening  may 
be  done  in  any  or  all  three  coordinate  directions.  Thus  there  exists  considerable 
generality  in  refinement  strategies  for  a  multiblock  grid.  However,  the  complete 
generality  of  a  true  local  grid  refinement  scheme  that  would  lead  to  a  more 
complicated  data  structure  does  not  exist. 

An  error  estimate  is  a  desirable,  though  often  unattainable,  component  in  an 
adaptive  grid  procedure.  Two  techniques  have  been  developed  to  deal  with  error 
estimation.  One  method  uses  the  grid  refinement  capability  previously  described. 
Under  certain  conditions,  the  error  in  a  numerical  solution  can  be  estimated  by 
comparing  fine  and  coarse  grid  solutions.  Another  method  of  error  estimation 
has  been  included.  The  method  is  based  on  the  premise  that  the  truncation  error 
in  a  difference  approximation  can  be  estimated  by  comparing  coarse  and  fine  grid 
difference  operators  evaluated  on  the  same  numerical  solution. 

Section  II  contains  a  discussion  of  grid  refinement  and  coarsening.  It  also 
includes  a  description  of  the  interpolation  schemes  that  are  used.  The  basic  idea 
and  the  governing  equations  that  lead  to  adaptive  grid  generation  are  presented 
in  Section  III.  Section  IV  discusses  the  development  of  several  grid  quality  mea¬ 
sures  and  the  truncation  error  estimation  techniques  for  use  in  the  construction 
of  adaptive  grid  generation.  A  discussion  of  the  results  of  static  adaptation  of 
several  grids  is  presented  in  Section  V.  Multiple  adaptation  performed  with  the 
adaptive  MISSE  Euler  flow  code  [3j  for  several  configurations  are  also  discussed 
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in  that  section.  Finally,  conclusions  and  recommendations  are  made  in  Section 
VI.  The  use  of  new  EAGLE  commands  and  namelists  along  with  several  exam¬ 
ples  are  given  in  Appendices  A  and  B.  The  procedure  for  the  incorporation  of 
adaptive  EAGLE  into  flow  codes  is  given  in  Appendix  C.  Other  improvements 
and  enhancements  that  have  been  made  on  the  flow  code  axe  given  in  Appendix 
D. 
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II.  GRID  REFINEMENT 

The  purpose  of  an  adaptive  grid  is  to  reduce  the  error  in  a  numerical  com¬ 
putation.  An  adaptive  grid  can  be  constructed  by  either  locally  refining  and/or 
coarsening  a  grid  or  by  moving  the  existing  grid  points.  A  local  grid  refine¬ 
ment  procedure  can  be  easily  incorporated  into  a  grid  or  flow  code  provided  the 
refinement  strategy  is  applied  to  complete  blocks. 

The  grid  refinement  procedure  involves  a  sequence  of  one-dimensional  in¬ 
terpolations.  Suppose  a  grid  is  defined  by  a  curvilinear  coordinate  system  with 
physical  variables  x,  y,  z  and  computational  variables  i  =  1,2,3.  By  splitt¬ 
ing  the  existing  physical  coordinate  values,  any  number  of  new  points  can  be 
inserted  by  evaluating  the  spline  function  at  new  values  for  the  curvilinear  coor¬ 
dinates.  The  new  points  will  lie  on  a  spline  curve.  It  should  be  noted  that  new 
boundary  grid  points  may  also  be  inserted  by  the  refinement  procedure.  These 
points  may  not  lie  on  the  actual  boundary  surface  that  was  used  to  construct 
the  original  grid.  However,  one  should  try  to  preserve  the  characteristics  of  the 
physical  boundary  as  far  as  possible.  For  this  reason  the  choice  of  the  spline 
function  is  important,  and  the  Akima  spline  has  been  used  in  this  work  to  avoid 
the  oscillations  that  may  occur  with  other  cubic  spline  methods.  In  the  process 
of  developing  an  efficient  grid  for  a  particular  problem,  it  may  be  decided  that 
some  of  the  grid  blocks  have  too  many  points  and  it  may  be  desirable  to  coarsen  a 
grid  in  certain  coordinate  directions.  This  can  also  be  done  as  a  sequence  of  one¬ 
dimensional  operations.  For  example,  along  the  grid  lines  in  a  single  coordinate 
direction,  every  other  grid  point  can  be  deleted.  This  can  then  be  repeated  along 
the  remaining  grid  lines  in  the  other  directions.  In  order  to  minimize  the  modi¬ 
fications  to  the  MISSE  flow  code,  the  cases  that  have  been  implemented  involve ' 


6 


only  halving  and  doubling  the  spacing  in  the  computational  region.  However, 
repeated  refinements  could  be  done  on  successive  restarts  of  a  numerical  solution. 

One  of  the  objectives  in  developing  the  grid  refinement  procedure  is  to  allow 
for  selective  refinement  and  coarsening  of  a  grid  as  a  solution  develops.  While 
this  has  not  been  automated,  the  user  has  the  option  of  examining  the  solution 
and  restarting  with  a  finer  or  coarser  grid  in  selected  blocks.  This,  in  particular, 
allows  for  an  initial  solution  to  be  started  on  a  coarse  grid  and  then  a  final  solution 
to  be  computed  on  a  refined  grid.  In  order  to  accomplish  this,  there  must  be 
the  capability  of  interpolating  between  coarse  and  fine  grids.  That  capability 
has  been  incorporated  into  the  flow  code.  Since  the  option  exists  for  refinement 
in  any  or  till  coordinated  directions,  the  interpolation  is  also  done  in  a  sequence 
of  one-dimensional  operations.  A  simple  linear  interpolation  in  computational 
space  is  used. 

Consider  a  function  <7*,  k  =  1,  •  •  • ,  kmax  defined  at  the  cell  centers  of  a  one¬ 
dimensional  grid.  On  a  grid  refinement,  the  new  cell  centers  are  located  at  points 
that  will  be  indexed  as  k  +  ~  and  k  +  | .  The  function  values  at  these  points  are 
computed  as 

3  1 

qk+\  =  4qk  + 

1  3 

9*+r  =  49*  +  4qk+1' 

With  a  cell-centered  solution  algorithm,  as  is  considered  here,  it  is  also  necessary 
to  use  interpolation  to  transfer  values  to  a  coarser  grid.  Reversing  the  above 
process,  if  a  coarse  grid  solution  value  <7*  is  needed,  it  can  be  computed  from  fine 
grid  solution  values  by  the  formula 

1, 

qk  =  I  +  7fe+i)- 

The  interpolation  procedure  is  repeated  in  each  coordinated  direction  in  which 
the  grid  has  been  refined  or  coarsened. 
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The  process  of  grid  refinement  and  coarsening  may  destroy  the  continuity 
of  grid  lines  at  block  boundaries.  Therefore,  a  similar  interpolation  scheme  has 
been  introduced  to  handle  the  block  boundary  interface. 


8 


III.  ADAPTIVE  GRID  GENERATION 
Variational  approach 

Among  the  various  adaptive  grid  generation  techniques,  the  minimization 
of  the  integral  of  some  grid  property  over  the  computational  domain  is  known 
as  the  variational  approach.  The  resulting  Euler  variational  equations  from  the 
calculus  of  variations  will  constitute  the  grid  generation  system.  The  choice  of 
what  property  is  to  be  minimized  depends  on  what  is  expected  from  the  grid.  For 
example,  Brackbill  and  Saltzman  [4]  constructed  this  adaptive  grid  by  minimizing 
a  weighted  combination  of  integrals  that  emphasize  smoothness,  orthogonality, 
and  concentration  of  the  grid.  A  similar  approach  considering  smoothness,  a 
measure  of  the  grid  cell  area,  and  the  orthogonality  of  the  grid  lines  can  be 
found  in  Roache  and  Steinberg  [5].  Several  other  grid  properties  that  might  be 
considered,  such  as  square  of  cell  volume  and  inverse  cell  volume,  can  be  found 
in  Warsi  and  Thompson  [6].  By  solving  a  large-scale  nonlinear  minimization 
problem  using  a  conjugate  gradient  method,  the  orthogonality  of  the  interior 
and  boundary  angles  has  been  controlled,  as  described  by  Castillo  [7]. 

Due  to  the  complexity  of  the  resulting  Euler  variational  equations,  they  are 
difficult  to  solve  and  solution  algorithms  may  not  converge.  Kennon  and  Du- 
likravich  [8]  and  Carcaillet,  et  al.  [9]  have  developed  algorithms  for  the  direct 
solution  of  the  variational  problem.  In  these  two  reports,  different  discrete  prob¬ 
lems  were  formulated  as  unconstrained  optimization  problems  and  then  solved 
by  a  conjugate  gradient  iterative  method.  A  survey  of  the  types  of  integrals  that 
may  be  included  in  a  variational  problem,  and  the  geometric  properties  that  each 
integral  imposes  upon  the  grid,  can  be  found  in  Soni  and  Mastin  [10]. 
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Control  function  approach 

The  control  function  approach  to  adaptation  is  developed  by  noting  that  the 
one-dimensional  form  of  the  elliptic  grid  generation  system  can  be  written  as 

x^  +  PXi  =  0,  (3-1) 

and  the  differentiated  form  of  the  equidistribution  principle,  Wx^  =  constant, 


Wx$£  +  W^x^  =  0, 


(3-2) 


where  P  is  the  function  to  control  the  coordinate  line  spacing,  and  W  is  the 
weight  function. 

From  (3.1)  and  (3.2),  the  control  function  can  be  defined  in  terms  of  the 
weight  function  and  its  derivative  as 


P  = 


Wz 

W  ' 


(3.3) 


This  equation  can  be  extended  in  a  general  three-dimensional  form  as 


P-^L 

w 


(3.4) 


This  approach  was  suggested  by  Anderson  [11,12]  and  has  been  applied  with 
success  for  two-dimensional  configurations  by  Johnson  and  Thompson  [13,14] 
and  for  three-dimensional  configurations  by  Kim  and  Thompson  [1]  and  by  Tu 
and  Thompson  [2]. 

The  complete  generalization  of  (3.4)  was  proposed  by  Eiseman  [15]  as 


p,  =  T9- 


9iJ  {Wi)v 
W,  ’ 


(3-5) 


where  Wi  is  the  weight  function  chosen  for  the  £*  direction.  This  definition  of  the 
control  functions  provides  a  convenient  means  to  specify  three  separate  control 
functions,  with  one  in  each  coordinate  direction. 


10 


In  order  to  preserve  the  geometrical  characteristics  of  the  existing  grid,  it 
is  practical  to  construct  the  control  functions  in  such  a  manner  that  the  control 
functions  defined  by  (3.5)  are  added  to  the  initial  set  of  control  functions  obtained 
from  the  geometry,  i.e., 

Pi*m9  +  Ci(Pi)wt  z  =  1,2,3,  (3.6) 

where 

(P,)g  :  control  function  based  on  geometry 
(P,)w  :  control  function  based  on  weight  function 
( Ci )  :  weight  coefficient  to  be  specified 

In  these  equations  the  weight  function  W  can  be  computed  by  different 
formulas  for  different  adaptive  mechanisms. 

Adaptation  to 

Variable  :  W  =  1+  |  V  \ 

Gradient  :  W  =  1+  |  W  |  (3-7) 

Curvature  :  W  =  (1  +  0  \  K  |)^1  +  a|  W  |2 

where  V  can  be  either  a  flow  solution  variable  or  a  grid  quality  measure.  Here 
a  are  in  the  interval  [0, 1]  and 

V2V 

I\  = - - - 5—  (3.8) 

(1  +  |  VV  |2)3/2 

is  the  curvature  of  the  variable  V. 

Using  these  definitions  of  the  control  functions,  the  elliptic  generation  system 
given  by 

i  =  1,2,3,  (3.9) 

or  transformed  to  the  computational  region  a s 

= °’ 

i=l 7=1  fc=I 


(3.10) 
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that  becomes  an  adaptive  grid  generation  system.  This  system  (3.10)  is  then 
solved  iteratively  in  adaptive  EAGLE  by  the  point  successive  overtaxation 
(SOR)  method  to  generate  the  adaptive  grid. 


12 


IV.  GRID  QUALITY  MEASURES  AND  SOLUTION  ERROR 

The  objective  of  this  part  of  the  investigation  is  to  develop  a  means  of 
evaluating  grids  through  the  computation  of  certain  grid  properties  that  are 
related  to  grid  quality  and  to  develop  techniques  for  estimating  the  truncation 
error. 

Following  Kerlick  and  Klopfer  [17]  and  Gatlin  [18],  the  grid  quality  measures 
are  taken  as  skew  Jingle,  aspect  ratio,  grid  Laplacian,  and  arc  length.  Techniques 
for  estimating  the  truncation  error  due  to  the  work  of  Mastin  [19]  are  also  in¬ 
cluded.  At  each  grid  point  in  a  general  three-dimensional  grid,  each  property 
can  have  three  values  associated  with  the  three  directions.  The  approach  taken 
under  this  investigation  is  to  treat  each  surface  of  constant  £*  separately  for  ease 
in  graphical  interpretation. 

Skew  angle 

The  minimum  skew  angle  between  intersecting  grid  lines  is  one  of  the  most 
important  measurable  grid  properties.  This  angle  can  be  expressed  in  terms  of 
the  covariant  metric  elements  as 


Since  g\%  =  g? i,  <713  =  <731,  and  <723  =  <732.  the  three  skew  angles  associated 
with  each  grid  point  in  a  three-dimensional  grid  are  #12,  #23.  and  #31.  Figure  1 
illustrates  the  skew  angle  at  a  grid  point  on  a  surface  of  constant  £3 . 

Aspect  ratio 

Since  aspect  ratio  is  the  ratio  of  the  lengths  of  the  sides  of  a  grid  cell,  it  can 
be  defined  in  two  different  ways.  For  example,  on  a  surface  of  constant  £3,  this 
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ratio  can  be  expressed  in  terms  of  metric  elements  gu  and  gjj  as 

ARij  =  (4.2a) 

V  9 a 

or 

ARji  =  ./^T  (4.26) 

V  0«» 

Large  changes  in  aspect  ratio  of  grids  from  one  part  of  the  field  to  another  may 
inhibit  the  convergence  of  viscous  flow  solutions  to  a  steady  state.  The  arc 
lengths  involved  in  the  geometric  representation  of  ARij  and  ARji  on  a  surface 
of  constant  £3  are  illustrated  in  Figure  2. 

Laplacian 

A  useful  measure  of  the  smoothness  of  a  grid  is  the  Laplacian  of  the  curvi¬ 
linear  system,  V2£‘,  i  =  1,2,3,  which  is  simply  the  rate  of  change  of  grid  point 
density  in  the  grid.  For  a  perfectly  uniform  grid,  the  grid  Laplacian  would  van¬ 
ish  everywhere,  but  exceedingly  large  values  may  arise  in  highly  stretched  grids. 
The  mathematical  representation  of  the  grid  Laplacian  is  defined  in  terms  of 
the  contravariant  metric  elements  g'* ,  the  contravariant  base  vectors  a1,  and  the 
position  vector  r  as 

VV  =  -£X>'^  •  »  -  1.2. 3-  (4-3) 

■=U=1 

Figures  3.1  and  3.2  illustrate  the  behavior  of  the  grid  Laplacian  on  a  surface  of 
constant  £3. 

Arc  length 

Another  important  measure  of  the  grid  quality  is  the  local  rate  at  which  grid 
spacing  changes.  On  a  coordinate  surface  of  constant  £3,  and  along  a  coordinate 
line  of  constant  £2 ,  the  grid  spacing  can  be  defined  els 


di  =  ((xi+i  -Xi)2  +  (yi+l  -yi)2  +(r<-H  -  St)2]*- 


(4.4a) 
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The  normalized  rate  at  which  grid  spacing  changes  (ARCL)  is  then 

(ARCL){  =  XTTTT^-  (4'44> 

2  (”*'  +  “i-1 ) 

Figure  4  illustrates  the  geometric  representation  of  ( ARCL)i . 

Error  estimation 

Finite  element  schemes  have  one  advantage  over  finite  difference  and  finite 
volume  schemes  in  that  rigorous  error  estimates  can  be  computed  in  many  cases. 
Although  truncation  error  estimates  can  be  derived  for  nearly  all  numerical  meth¬ 
ods,  they  frequently  involve  higher  order  derivatives  of  the  solution  that  cannot 
be  estimated  with  any  accuracy  using  the  numerical  solution.  A  different  ap¬ 
proach  is  taken  here.  Based  on  truncation  error  estimates  that  are  computed  in 
the  multigrid  methods,  the  truncation  error  is  estimated  by  applying  a  coarse 
and  a  fine  grid  difference  approximation  to  the  numerical  solution.  Another, 
more  rigorous,  approach  to  estimating  error  is  to  compare  coarse  and  fine  grid 
solutions.  Such  a  procedure  has  long  been  used  in  the  numerical  solution  of  ordi¬ 
nary  differential  equations  and  is  referred  to  as  Richardson  extrapolation.  This 
approach  can  now  be  carried  out  using  the  grid  refinement  procedures,  which 
have  been  included  in  the  MISSE  flow  code.  Both  of  these  techniques  should  be 
considered  as  heuristic  error  estimates,  which  are  intended  only  to  give  tin  order 
of  magnitude  estimates  of  the  truncation  error  and  the  solution  error. 

The  derivation  that  follows  can  be  applied  to  the  Euler  equations  for  com¬ 
pressible  flow  and  other  systems  of  conservation  laws  of  the  form 

+  ft  +  9y  +  =  0-  (4.5a) 

The  transformation  of  this  system  to  an  arbitrary  curvilinear  coordinate  system 
is 


Ut  +  (F?i  +  +  H^a)  =  0, 


(4.56) 
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where  y/g  is  the  Jacobian  of  the  transformation  and 

F  =  y/g((lf  +  (U  +  (',h) 

G=vn(if+e,s+e,k) 

H  =  +  &  +  £'■)• 

Let  h  be  the  spacing  of  the  fine  grid  and  nh  be  the  spacing  of  the  coarse  grid. 
Let  Ch  be  the  difference  approximation  operator  of  order  p  on  the  fine  grid  and 
Cnh  be  the  same  difference  approximation  operator  on  the  coarse  grid.  Then  the 
finite  difference  approximation  of  the  PDE  can  be  represented  on  the  fine  grid  as 

«,  +  /,  +  g,  +  h,  =  Ch(F ;  G,  H)  +  T(h)>  (4.5  c) 

and  on  the  coarse  grid  as 

+  /«  +  3,  +  ht  =  Cnh(F,G,H)  +  T{nhY ,  (4.5d) 

where  n  is  an  integer.  From  (4.5c)  and  (4.5d),  the  estimate  of  the  truncation 
error  on  the  fine  grid  can  be  computed  as 

W  =  " f  ■  (4-6) 

Richardson  extrapolation  can  be  used  to  compute  the  error  in  the  numerical 
solution.  Even  though  numerical  solutions  must  be  computed  on  both  fine  and 
coarse  grids,  the  error  estimates  that  result  do  not  have  the  large  peaks  at  so¬ 
lution  singularities  that  can  be  encountered  with  the  truncation  error  described 
above.  Thus  the  solution  error  estimates  may  sometimes  be  more  useful  in  the 
construction  of  adaptive  grids. 

Assume  that  there  are  two  numerical  solutions  of  order  p  accuracy  for  (4.56) 
that  have  been  computed  on  a  fine  grid  and  on  a  coarse  grid,  with  grid  spacing 
h  and  nh,  respectively,  in  each  coordinate  direction.  Assuming  that  the  same 
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pth  order  method  is  used  in  both  cases,  the  relation  between  the  two  numerical 
solutions  and  the  actual  solution  u  of  the  PDG  can  be  established  as 


u  =  Uh  +  R{h)p 


(4-7  a) 


and 

u  =  Unh  +  R(nh)p.  (4.7  b) 

Prom  these  equations,  an  extrapolated  value  of  u  can  be  computed  as 

npUk  -  Unh 


u  = 


(4.7c) 


(nP-l) 

Thus  the  estimate  of  the  error  in  the  numerical  solution  computed  on  the  h  grid 


is: 


u-Uk  = 


Uk  ~  Unk 

(up  -  1)  ‘ 


(4.7d) 
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V.  RESULTS  AND  DISCUSSION 

The  adaptive  grid  generation  system  based  on  the  control  function  approach, 
as  described  in  the  previous  sections,  has  been  used  to  generate  static  and  mul¬ 
tiple  adaptive  grids  for  several  geometries.  Results  are  presented  in  this  section. 
The  static  adaptive  grids  were  obtained  by  adapting  the  initial  grids  to  either 
grid  quality  measure  variables  or  to  existing  flow  solution  variables.  The  mul¬ 
tiple  adaptive  procedure  was  tested  on  several  different  configurations  with  the 
adaptive  MISSE  Euler  flow  code  for  transonic  and  supersonic  flow  cases. 

Static  adaptation  results 

Some  examples  of  the  grid  quality  adaptation  are  shown  in  Figures  5  through 
13,  for  adaptation  to  various  quality  measures.  Contour  plots  of  the  quality 
measures  are  given  for  comparison  with  those  for  the  initial  elliptic  grid  shown 
in  Figure  5.  (In  all  these  plots,  the  colors  range  from  blue  through  green  to  red 
as  the  value  of  the  quantity  represented  increases.)  In  these  contour  plots  the 
Laplacian  is  a  smoothness  measure,  with  blue  being  the  smoothest.  In  the  skew 
angle  plots,  red  is  the  most  orthogonal.  In  the  aspect  ratio  plots,  green  is  the 
closest  to  unity,  while  blue  and  red  indicate  large  ratios  in  the  two  directions. 
Approximately  the  same  number  of  total  adaptive  iterations  were  run  in  each 
case.  The  value  given  for  the  variable  ITMAX  on  each  figure  is  the  number  of 
inner  iterations  between  reevaluations  of  the  quality  and  application  of  further 
adaptation.  The  value  for  variable  ITMAXA  is  the  total  number  of  adaptations. 
The  value  for  AFIXP  =  NO  means  the  control  functions  Eire  updated  based  on 
the  geometry  of  the  previous  grid  rather  than  the  initial  grid,  at  each  adaptation. 
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Static  adaptation  to  quality  measures 

Comparison  of  Figure  6  with  Figure  5  shows  that  adaptation  to  the  skewness 
is  effective  in  reducing  the  skewness  in  one  region,  while  increasing  the  skewness 
in  other  regions  of  the  grid.  A  small  improvement  in  aspect  ratio  occurs,  but  the 
smoothness  of  the  grid  is  decreased.  Figure  7  shows  the  difference  of  the  average 
skew  angle  between  the  initial  and  adaptive  grids. 

Comparison  of  Figure  8  with  Figure  5  indicates  that  adaptation  to  aspect 
ratio  does  improve  both  aspect  ratio  and  smoothness  of  the  grid;  the  skewness 
is  increased,  however.  Figure  9  shows  the  large  change  in  average  skew  angle 
for  this  case.  Comparison  of  Figure  10  with  Figure  5  shows  that  adaptation  to 
smoothness  improves  the  skewness  and  aspect  ratio  of  the  grid  effectively,  but 
the  adaptive  grid  is  not  as  smooth  as  the  initial  elliptic  grid.  The  average  skew 
angle  of  the  initial  and  adaptive  grids  for  this  case  are  plotted  in  Figure  11. 

Figures  7,  12,  and  13  show  the  beneficial  effect  of  including  adaptation  to 
aspect  ratio,  arc  length,  and  smoothness  with  adaptation  to  skewness:  the  skew¬ 
ness  is  reduced  more  by  the  combination  than  with  skewness  adaptation  alone. 
A  little  improvement  occurs  in  aspect  ratio;  the  smoothness  of  the  grid  does, 
however,  decrease. 

Results  from  these  examples  show  that  the  adaptation  to  the  combination  of 
all  grid  quality  measures  or  to  each  individually  can  improve  some  grid  properties 
while  damaging  others.  For  example,  the  adaptation  to  the  Laplacian  of  this 
particular  grid  can  reduce  the  skewness,  but  the  resulting  adaptive  grid  is  not  els 
smooth  as  the  initial  grid.  The  choice  of  the  adaptive  variable  for  the  adaptation 
very  much  depends  on  what  property  of  the  grid  needs  to  be  improved  and  the 
configurations  of  the  grids. 

Static  adaptation  to  flow  solution 

Another  new  feature  of  the  adaptive  EAGLE  grid  code,  in  addition  to  the 
static  adaptation  to  the  quality  measure,  is  the  adaptation  to  an  existing  flow 


solution.  Examples  demonstrating  this  new  feature  are  shown  in  Figures  14 
through  16.  A  value  given  for  ITMAXA  =  3  on  each  plot  indicates  the  total 
number  of  adaptations,  and  AFIXP  =  NO  indicates  the  control  functions  were 
computed  based  on  the  geometry  of  the  previous  grid.  The  density  variable  from 
the  Euler  solution  (MISSE  Euler  code)  with  a  Mach  number  of  2  on  a  double 
wedge  (81  by  20)  grid  is  considered  in  this  case.  Figure  14  shows  the  color 
contour  plot  of  the  density  obtained  from  a  200-time-step  solution  on  the  initial 
grid.  This  grid  was  then  adapted  statically  to  the  gradient  of  the  density  in  one 
direction  with  a  weight  coefficient  Ci  of  1.5  and  to  the  curvature  of  the  density 
in  the  other  direction  with  a  weight  coefficient  C2  of  0.5,  a  of  1,  and  j3  of  0.5, 
(AWT=GRAD,CURV,  and  RHO=l,l).  The  color  plots  of  the  control  functions 
and  weight  functions  in  each  direction  are  also  shown  in  this  figure. 

Figure  15  shows  another  case  of  the  grid  being  adapted  statically  to  the 
curvature  of  the  density  in  one  direction  with  a  C\  of  0.5,  a  of  1,  and  0  of  0.5, 
and  to  the  gradient  of  the  density  in  the  other  direction  with  a  C2  of  1.5.  Figure 
16  shows  the  initial  grid  and  several  adaptive  grids  which  result  from  different 
adaptive  mechanisms  and  parameters. 

There  are  several  significant  advantages  from  these  new  features  of  the  adap¬ 
tive  EAGLE  grid  code  that  cannot  be  found  in  the  earlier  form  [1,2].  For  example, 
grid  quality  measures  along  with  control  and  weight  functions  in  each  direction 
of  the  resulting  generated  grid  can  be  output  for  graphical  representation.  This 
makes  it  easier  for  EAGLE  users  to  properly  examine  each  property  of  the  re¬ 
sulting  generated  grid.  Moreover,  these  grid  quality  measures  can  be  used  in 
addition  to  the  flow  solution  in  the  weighted  average  during  the  multiple  adap¬ 
tation  process.  The  capability  of  static  adaptation  to  an  existing  flow  solution  of 
the  adaptive  EAGLE  grid  code  is  sometimes  very  helpful  in  the  multiple  adap¬ 
tive  process.  Since  it  is  not  clear  what  weight  functions  and  weight  coefficients 
or  what  flow  variable  or  how  many  adaptive  iterations  will  produce  a  reasonable 
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grid  for  different  type  of  flows,  the  users  can  use  the  static  adaptive  mode  of  the 
EAGLE  grid  code  to  determine  these  parameters.  The  determination  of  these 
parameters  prior  to  submitting  a  multiple  adaptive  run  of  the  adaptive  flow  code 
can  save  the  users  a  considerable  amount  of  time. 

Multiple  adaptation  results 

Results  of  multiple  adaptation  performed  with  the  adaptive  MISSE  Euler 
flow  code  are  shown  in  Figures  17  through  34.  In  ail  these  plots,  NIT  is  the 
total  number  of  time  steps,  INT  indicates  the  number  of  time  steps  at  which 
the  first  adaptation  is  performed,  NCL  is  the  number  of  time  steps  between 
each  adaptation,  and  MAXINT  indicates  the  number  of  time  steps  at  which 
the  last  adaptation  is  performed.  Values  of  weight  functions  AWTi,  AWT2; 
weight  coefficients  C\,C2i  adaptive  variables  density  RHOj,  RHO2;  pressure 
PRESi,  PRES2  are  given  for  and  £2  directions,  respectively  For  example, 
AWT=GRAD,  CURV,  Cx  =  .5,  C2  =  .3,  RHO=1,0,  PRES=0,1,  and  a  =  1,0  = 
1  can  be  interpreted  as  the  adaptation  to  the  density  gradient  in  f1  direction 
with  a  Ci  of  .5  and  to  the  curvature  of  the  pressure  in  £2  direction  with  a  C2 
of  .3  and  with  coefficients  of  gradient  and  curvature  being  a  =  1  and  3=1. 
respectively. 

Multiple  adaptation  for  the  airfoil  RAE  2822 

The  first  case  considered  here  is  a  transonic  airfoil  (RAE  2822)  with  3.22 
degree  angle  of  attack  and  Mach  number  of  0.728.  Figure  17(c)  shows  the  density 
contours  obtained  from  500  time  steps  on  the  initial  grid.  (81  by  41)  while  the 
density  contours  obtained  from  the  same  number  of  time  steps  on  the  adaptive 
grid  are  in  17(d).  The  adaptation  is  to  the  gradient  of  the  combination  density 
and  pressure  in  both  directions  ( AWT=GRAD,GRAD,  RHO=l,l,PRES=l,l) 
with  a  Ci  of  .75,  and  a  C2  of  .2.  The  total  number  of  adaptations  is  5. 

Figure  18  shows  the  pressure  contours  on  the  initial  and  the  adaptive  grids. 
The  pressure  coefficients  of  the  upper  and  lower  surfaces  of  this  airfoil  are  plotted 
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in  Figure  19  along  with  the  experimental  data.  From  this  figure  the  adaptive 
grid  gives  a  better  prediction  of  the  location  of  the  shock  and  a  much  steeper 
shock  than  the  initial  grid.  A  record  of  the  CPU  time  on  an  IRIS  4D/440VGX 
machine  shows  that  the  total  CPU  time  for  the  initial  grid  without  adaptation 
was  1499.93  CPU  seconds  and  for  the  adaptive  grid  was  1518.63  CPU  seconds, 
an  1.2  percent  increase. 

Multiple  adaptation  for  a  double  wedge 

The  major  goal  of  multiple  adaptation  is  to  move  the  grid  points  to  capture 
major  features  of  the  flow  field  as  the  flow  develops.  In  supersonic  flow,  these 
features  are  shocks  and  expansions.  Flow  over  a  double  wedge  and  into  a  wind 
tunnel  are  considered  in  the  present  investigation. 

Results  obtained  from  a  supersonic  flow  at  a  Mach  number  of  2  over  fine 
(121  by  41)  and  coarse  (81  by  31)  double  wedge  grids  are  shown  in  Figures  20 
through  32.  Figure  20  shows  the  density  contours  obtained  from  300  time  steps 
on  the  initial  and  adaptive  grids  (121  by  41).  The  grid  was  adapted  to  the  density 
gradient  in  the  flow  direction  (RHO  =  1,0)  with  a  C\  of  .7  and  to  the  pressure 
gradient  in  the  normal  direction  with  a  Ci  of  .5.  A  total  of  four  adaptations  was 
used  for  this  case,  and  the  variable  AFIXP  =  NO. 

Figure  21  shows  the  pressure  contours  on  the  initial  and  adaptive  grids 
(121  by  41).  Figure  22  shows  the  pressure  coefficients  of  the  lower  wall  and  a 
convergence  history  of  the  two  solutions  are  shown  in  Figure  23.  In  Figure  23, 
the  high  peaks  at  each  adaptation  are  due  to  the  use  of  the  previous  solutions 
on  the  new  adapted  grid.  From  these  figures,  clearly  the  adaptive  grid  gives  a 
much  better  representation  of  the  shock  regions  as  well  as  the  expansion  regions. 
Shocks  are  much  sharper  for  the  solution  obtained  on  the  adaptive  grid.  A 
record  of  the  CPU  time  on  an  IRIS  4D/440VGX  machine  shows  that  the  total 
CPU  time  for  the  initial  grid  (121  by  41)  without  adaptation  was  1481.51  CPU 
seconds  and  for  the  adaptive  grid  (121  by  41)  was  1599.02  CPU  seconds,  an  8 
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percent  increase.  With  this  time  record,  it  is  worthy  to  use  multiple  adaptation 
to  achieve  a  better  representation  of  the  solution.  Contour  plots  of  the  density 
for  the  initial  fine  grid  (121  by  41),  initial  coarse  grid  (81  by  31),  and  adaptive 
grid  (81  by  31)  are  shown  in  Figure  24.  In  this  figure,  solutions  on  the  right 
panel  were  obtained  from  the  corresponding  grids  on  the  left  panel.  The  coaxse 
grid  was  adapted  to  the  combination  of  density  and  pressure  in  f 1  direction  with 
weight  coefficient  C\  of  .5  and  to  the  gradient  of  this  combination  in  f2  direction 
with  weight  coefficient  C2  of  .5,  ( AWT = VAR, GRAD,  RH0=1,1,  PRES=1,1). 
The  total  number  of  adaptations  is  four  with  AFIXP=NO.  Figure  25  shows  the 
pressure  contours  obtained  on  each  grid.  Pressure  coefficients  of  the  lower  walls 
obtained  from  the  fine  and  coarse  grids  are  plotted  in  Figure  26.  Figure  27  shows 
the  pressure  coefficients  of  the  lower  walls  of  the  fine  and  adaptive  coarse  grids. 

Different  adaptive  mechanisms  were  applied  to  the  coarse  grid  in  the  multiple 
adaptation  process  are  shown  in  Figures  28  through  33.  Figure  28(b)  shows  the 
density  contours  obtained  on  the  adaptive  grid  of  Figure  28(a).  The  initial  grid 
was  adapted  to  the  curvature  of  the  combination  of  density  and  pressure  in 
both  directions  ( AWT=CURV,CURV,  RHO=l,l,  PRES=1,1).  Total  number  of 
adaptations  is  four  with  aCi  of  .7  and  a  C2  of  .7.  The  coefficients  of  the  gradient 
and  curvature  are  a  =  1  and  0  =  .5,  respectively,  and  AFIXP=YES.  Figure  30 
shows  the  pressure  coefficient  of  the  lower  wall  for  this  case. 

Figure  28(d)  shows  the  density  contours  obtained  on  the  adaptive  grid  of 
Figure  28(c).  The  adaptive  mechanism  for  this  case  was  pressure  gradient  in 
both  directions  with  a  C]  of  .7,  a  C2  of  .7,  and  total  number  of  adaptations  is 
four  ( AWT=GRAD.GRAD,  RHO=0.0,  PRES=1,1).  The  pressure  coefficient  on 
the  lower  wall  obtained  from  this  adaptation  is  shown  in  Figure  31. 

The  initial  grid,  adapted  to  the  gradient  of  the  combination  of  density  and 
pressure  in  £2  direction  only,  is  shown  in  Figure  28(e).  Total  number  of  adap¬ 
tations  is  five  with  a  C\  of  0.,  a  C2  of  .9,  and  AFIXP  =  NO.  Density  contours 
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obtained  from  this  adaptive  grid  are  shown  in  Figure  28(f),  and  the  pressure 
coefficient  on  the  lower  wall  is  plotted  in  Figure  32.  Adaptive  grids  and  pressure 
contours  of  these  adaptive  mechanisms  are  shown  in  the  left  and  right  panel  of 
Figure  29,  respectively. 

From  these  figures,  the  representation  of  the  shocks  on  the  adapted  coarse 
grid  is  much  sharper  and  closer  to  the  fine  grid  solution  than  the  nonadap- 
tive  coarse  grid.  The  total  CPU  time  for  obtaining  300  time  steps  solution  for 
the  adaptive  grid  was  approximately  800  seconds  for  each  adaptive  mechanism, 
nearly  a  50  percent  time  saving  compared  to  that  of  the  fine  grid. 

The  adaptation  to  the  combination  of  density  and  pressure  in  f1  direction 
and  to  the  gradient  of  this  combination  in  £2  direction  of  Figure  24  gives  a 
smoother  pressure  coefficient  behind  the  shock  than  the  adaptation  to  the  gradi¬ 
ent  of  pressure  alone  of  Figures  28(c)  and  (d).  The  adaptation  to  the  curvature  of 
Figures  28(a)  and  (b)  gives  a  better  result,  however,  with  a  little  over  prediction 
of  the  pressure  coefficient  right  behind  the  shock.  The  adaptation  to  the  gradient 
of  the  combination  of  the  density  and  pressure  in  <f2  direction  only  of  Figures 
28(e)  and  (f)  gives  the  closest  solution  to  the  fine  grid  solution.  The  convergence 
history  of  solutions  obtained  from  this  adaptive  grid  and  non-adaptive  fine  and 
coarse  grids  are  plotted  in  Figure  33. 

From  these  results,  clearly  multiple  adaptive  grids  produce  a  better  rep¬ 
resentation  of  the  shock  regions  as  well  as  the  expansion  regions  than  that  of 
the  same  nonadaptive  grid.  Among  these  adaptive  mechanisms,  the  use  of  the 
weighted  average  of  weight  functions  computed  from  several  flow  variables  gives 
better  results  than  the  use  of  a  single  variable.  Another  advantage  that  should  be 
mentioned  here  is  the  controlling  of  the  direction  in  which  adaptation  is  applied. 
As  shown  above,  the  adaptation  in  only  one  direction  (£2)  gives  the  closest  solu¬ 
tion  to  the  fine  grid  solution.  Moreover,  the  grid  in  this  adaptive  mechanism  is 
not  disturbed  as  much  as  the  adaptation  in  both  directions.  The  minimum  skew 
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angle  for  this  case  is  higher  compared  to  those  of  adaptation  in  both  directions. 
Of  course,  this  is  true  only  for  a  certain  number  of  adaptations  and  a  particular 
value  of  weight  coefficients. 

Multiple  adaptation  for  a  wind  tunnel 

Results  from  the  supersonic  flow  at  a  Mach  number  of  2  in  a  wind  tunnel 
are  shown  in  Figures  34  through  35.  These  results  were  also  obtained  in  300  time 
steps.  Figure  34(a)  is  the  initial  grid,  (c)  is  the  adaptive  grid  adapted  to  the  error 
estimation  in  both  directions,  and  (e)  is  the  adaptive  grid  adapted  to  gradient 
of  the  combination  of  density  and  pressure  in  both  directions.  The  number  of 
adaptations  is  five  for  both  cases,  with  a  Ci  of  1  and  a  C2  of  1  for  adaptation  to 
error  estimation,  and  a  C\  of  0.6  and  a  Ci  of  0.55  for  the  adaptation  to  gradient 
of  the  combination.  Figure  34(b),  (d),  and  (f)  are  the  pressure  contours  on  initial 
and  adaptive  grids.  Figure  35  shows  the  contour  plot  of  the  error  estimation  on 
each  grid.  Shocks  are  much  sharper  for  solutions  obtained  on  the  adaptive  grids 
than  the  nonadaptive  grids  for  this  configuration  in  supersonic  flow  as  well. 

Results  from  these  examples  show  that  multiple  adaptive  grids  performed 
well  in  capturing  major  features  of  the  flow  field  in  supersonic  flow  for  these 
particular  configurations.  The  adaptations  to  the  combination  of  the  grid  quality 
measures  such  as  skewness  of  the  grid  and  the  flow  solution  for  these  particular 
grids  not  only  make  the  grid  more  skewed  but  also  resulted  in  poor  resolution 
of  the  major  features  of  the  flow  field.  On  the  other  hand,  the  adaptation  to 
the  error  estimation  and  the  use  of  the  weighted  average  in  weight  functions 
computed  from  several  flow  variables  does,  in  fact,  improve  the  solutions. 

The  computation  of  the  weight  functions  and  the  choice  of  the  adaptive 
solution  variable  are  independent  from  one  direction  to  another,  which  enable  the 
users  to  have  more  freedom  in  choosing  suitable  adaptive  mechanism  for  each  kind 
of  flow.  For  example,  in  the  case  of  boundary  layers  and  shocks  occurring  in  the 
same  flow  field,  the  users  may  choose  to  adapt  'the  grid  to  the  velocity  magnitude 
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gradient  in  the  normal  direction  to  capture  the  boundary  layer  regions  and  to 
the  pressure  gradient  in  the  flow  direction  to  capture  the  shocks. 

Grid  refinement 

Block  grid  refinement  can  be  used  to  improve  the  quality  of  numerical  solu¬ 
tions.  In  the  following  examples  a  single  block  grid  is  divided  into  subblocks  to 
demonstrate  this  concept.  It  should  be  noted  that  the  subdivision  and  refinement 
strategy  is  done  manually  for  these  examples.  This  could  be  done  automatically 
for  simple  configurations  like  the  examples  considered  here.  However,  additional 
research  is  still  needed  before  the  process  can  be  automated  for  more  general 
problems. 

The  first  example  is  the  solution  about  an  airfoil.  The  solution  is  for  a  Mach 
number  of  0.8  and  a  1.25  degree  angle  of  attack.  The  regions  of  interest  for 
this  problem  are  the  stagnation  point  at  the  leading  edge,  a  strong  shock  on  the 
upper  surface,  and  a  weak  shock  on  the  lower  surface.  A  solution  on  a  81  by 
41  grid  is  plotted  in  Figure  36.  The  truncation  error  estimate  for  the  solution 
vector  is  plotted  in  Figure  37.  The  error  estimate  is  effected  in  this  example  by 
difficulties  near  the  airfoil  surface  caused  by  characteristic  boundary  conditions. 
The  grid  is  then  subblocked  and  refined  in  both  directions  in  the  region  around 
the  shock  on  the  upper  surface  and  near  a  possible  shock  on  the  lower  surface. 
The  solution  computed  on  this  locally  refined  grid  shows  a  marked  difference  in 
the  location  and  strengths  of  the  shocks,  as  seen  in  Figure  38.  This  solution  is 
compared  with  the  refined  single  block  161  by  81  solution  appearing  in  Figure 
39.  The  two  solutions  are  nearly  identical,  indicating  that  the  local  block  refined 
solution  is  as  accurate  as  the  fine  grid  solution.  There  is  no  noticeable  effect  of  the 
interpolation  at  the  interfaces.  The  coarse  and  fine  grid  solutions  in  Figures  36 
and  39  were  used  to  compute  and  estimate  the  solution  error,  using  Richardson 
extrapolation.  That  estimate  is  plotted  in  Figure  40. 


26 


The  next  application  of  block  grid  refinement  is  for  flow  about  a  three- 
dimensional  cylindrical  store.  The  solution  is  computed  with  a  Mach  number  of 
1.2  and  a  zero  angle  of  attack.  The  first  solution  is  plotted  in  Figure  41.  The 
shock  structure,  which  should  form  near  the  front  and  rear  of  the  body,  is  hardly 
distinguishable  due  to  the  coarseness  of  the  grid.  These  features  can  be  brought 
out  by  refining  the  grid,  as  seen  in  Figure  42.  The  truncation  error  computed 
from  the  coarse  grid  solution  is  also  plotted  in  Figure  43.  The  largest  values  occur 
near  the  front  and  rear  stagnation  points  and  also  along  the  axis  of  symmetry 
where  the  Jacobian  vanishes. 
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VI.  CONCLUSIONS 

The  widely  used  EAGLE  grid  generation  system  has  been  extended  and 
enhanced  so  that  it  can  be  readily  coupled  with  existing  PDE  solvers  that  operate 
on  structured  grids  to  provide  a  flexible  adaptive  grid  capability.  The  adaptive 
EAGLE  grid  code  now  can  be  used  for  generating  not  only  algebraic  grids  and 
elliptic  grids  but  static  adaptive  grids  as  well.  In  the  static  adaptation,  the  grid 
can  be  adapted  to  an  existing  PDE  solution  or  to  grid  quality  measures  or  to 
a  combination  of  both.  The  test  cases  show  that  some  grid  properties  can  be 
improved  by  the  static  adaptation  to  grid  quality  measures. 

In  this  study,  the  weight  functions  can  be  formulated  as  weighted  average 
of  weight  functions  from  several  flow  variables  or  several  quality  measures  or 
the  combination  of  both.  Different  weight  functions  and  adaptive  variables  can 
be  applied  in  each  direction.  These  operations  are  controlled  through  the  input 
parameters  in  static  as  well  as  multiple  adaptation  mode. 

There  are  several  successful  incorporations  of  the  adaptive  EAGLE  subrou¬ 
tines  into  flow  codes,  including  the  MISSE  Euler  solver.  Several  configurations 
are  considered  for  each  of  these  adaptive  flow  codes  for  the  investigation  of  the 
new  weight  functions  formulations  and  grid  quality  measures  in  the  multiple 
adaptation.  Results  obtained  from  the  adaptive  MISSE  Euler  flow  code  show 
considerable  success  as  measured  by  improvements  in  shock  resolution  on  coarse 
grids  in  the  compressible  flows. 

Several  enhancements  to  the  MISSE  code  itself  will  aid  in  the  estimation  and 
control  of  error.  These  include  the  capability  of  block  refinement  and  coarsening 
and  the  computation  of  a  truncation  error  estimate. 
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APPENDIX  A 

EAGLE  ADAPTIVE  COMMANDS 

NAME 

RESTART  generates  grid  from  restart  file. 

SYNOPSIS 

$  ’RESTART’,  FILNAM  = _ _  VARIN  = - $ 

DESCRIPTION 

The  function  of  this  command  is  to  read  in  the  restart  file  from  the  previous 
run.  This  file  is  unformatted  and  its  name  has  been  set  to  ’rsfile’. 

PARAMETER 

FILNAM  =  ’rsfile’ 

VARIN  =  ’NO’  indicates  that  during  the  previous  run  the  adaptive  variables 
array  has  not  been  saved  on  the  restart  file. 

VARIN  =  ’YES’  indicates  that  during  the  previous  run  the  adaptive  vari¬ 
ables  array  has  been  saved  on  the  restart  file,  and  adaptation  is  to  be  done  on 
the  current  run  based  on  this  adaptive  variables  array. 

EXAMPLE 

$  ’BLOCK’  ,  SIZE=35,21,1$ 

$  ’BLOCK’  ,  SIZE=80,41,1$ 

$  ’RESTART’  ,  FILNAM  =  ’rsfile’  ,  VARIN  =  ’NO’  $ 

This  grid  runstream  begins  with  the  command  BLOCK’,  which  indicates 
that  the  previous  runstream  generated  a  two-block  grid.  The  first  block  has 
dimensions  35  x  21  x  l.and  the  second  block  has  dimension  80  x  41  x  1.  The 
order  and  the  size  of  each  block  in  a  restart  runstream  must  be  consistent  with 
that  of  the  previous  one.  The  adaptive  variables  array  is  not  contained  in  this 
’rsfile’,  as  indicated  by  VARIN  =  ’NO’. 
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See  also: 

RESTART  =  ’YES’ 

Bugs: 

Errors  may  occur  if  the  namelists  VAROUT  and  VARIN  are  not  specified 
properly. 

Notes: 

The  word  RESTART  here  is  used  in  two  cases:  the  namelist  name  RESTART 
and  the  command  RESTART.  The  same  word  has  two  different  functions. 

NAME 

VFILE  reads  in  the  flow  solution  variables  file. 

% 

SYNOPSIS 

$  ’VFILE’,  FILNAM  = _ _  FORM  = _ $ 

DESCRIPTION 

The  function  of  this  command  is  to  read  in  the  flow  solution  variables  file 
to  perform  static  adaptation. 

PARAMETER 

FILNAM  is  the  name  of  the  solution  file. 

FORM  =  ’LIST’  indicates  the  form  of  this  file  is  formatted,  and  it  has  the 
same  format  as  the  Q  file  of  PL0T3D. 

EXAMPLE 

$  ’VFILE’  ,  FILNAM  =  ’soln.fmt’  ,  FORM  =  ’LIST’  $ 

Bugs: 

To  avoid  errors  in  this  case,  this  line  should  always  be  right  before  the  CUT 
command,  or  otherwise  before  the  command  END  of  input. 
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APPENDIX  B 

EAGLE  ADAPTIVE  NAMELISTS 

PARAMETER 

ITMAXA  is  the  number  of  adaptive  iterations. 

PARAMETER 

ADAPT  =  ’NO’  or  ’NONE’  sure  default  values;  elliptic  grid  is  produced. 
ADAPT  =  ’YES’  produces  the  adaptive  grid. 

PARAMETER 

AWT  =  ’VAR’,  ’VAR’,  ’VAR’  is  the  adaptation  to  variable  with  the  weight 
function  as 

W  =  1+|  V  |, 

where  V  is  either  a  flow  variable  or  a  quality  measure  variable. 

AWT  ==  ’GRAD’,  ’GRAD’,  ’GRAD’  is  the  adaptation  to  the  gradient  of  the 
variable  with  the  weight  function  as 

W  =  1+  |  VV  |, 

where  V  is  either  a  flow  variable  or  a  quality  measure  variable. 

AWT  =  ’CURV’,  ’CURV’,  ’CURV’  is  the  adaptation  to  the  curvature  of  the 
variable  with  the  weight  function  as 

W  =  (1  +  0  |  K  \)yj l+a|  W  |2, 

where 

K  _  v2v 

(l  +  |VV(2)3/2- 

V  is  either  a  flow  variable  or  a  quality  measure  variable. 


34 


Notes: 

One  may  specify  AWT  =  'VAR’,  ’GRAD’,  ’CURV1,  which  means  the  grid  is 
adapted  to  the  variable  in  the  1  direction,  to  the  gradient  in  the  2  direction  and 
to  the  curvature  in  the  3  direction.  Any  combination  of  this  is  also  valid. 

PARAMETER 

CW:  Weight  coefficients,  default  values  are  1.0,  1.0,  1.0. 

PARAMETER 

ALPHA:  Coefficient  of  the  gradient,  in  the  range  form  0  to  1,  default  to  1.0. 

PARAMETER 

BETA:  Curvature  coefficient,  in  the  range  from  0  to  1,  default  to  1.0. 
Notes: 

The  following  set  of  parameters,  ASKEW  through  VORR,  represent  the 
geometric  and  solution  variables  for  the  adaptive  process.  The  default  value  is 
0.0.  A  value  of  1.0  indicates  the  variable  is  used  in  the  calculation  of  the  weight 
function. 

PARAMETER 

ASKEW:  Skew  angle. 

PARAMETER 

AASPE:  Aspect  ratio. 

PARAMETER 

AARCL:  Arc  length. 

PARAMETER 


APLAC:  Laplacian. 
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Notes: 

If  ail  of  these  namelists  are  not  0.0,  then  the  weight  function  is  computed  as 
a  weighted  average  of  the  individual  weight  functions. 

PARAMETER 

RHO:  Density. 

PARAMETER 

RHOU:  X- momentum. 

PARAMETER 

RHOV:  Y-momentum. 

PARAMETER 

RHOW:  Z-momentum. 

PARAMETER 

RHOE:  Energy. 

PARAMETER 

VOMA:  Velocity  magnitude. 

PARAMETER 

VORR:  Vorticity  magnitude. 

PARAMETER 

VARIN:  Defaulted  value  is  ’NO’,  VARIN  =  ’YES’  indicates  the  restart  file 
’rsfile'  contains  adaptive  variables  array  from  the  previous  run. 

PARAMETER 

VAROUT:  Defaulted  value  is  ’NO’,  VAROUT  =  ’YES’  indicates  that  in  the 
current  run  adaptive  variables  array  will  be  saved  on  the  restart  file  ’rsfile’. 
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PARAMETER 

RESTART:  Defaulted  value  is  ’NO’,  RESTART  =  ’YES’  means  a  restart 
file,  namely  ’rsfile’,  will  be  generated  at  the  end  of  the  current  run. 

PARAMETER 

AFIXP:  Defaulted  value  is  ’YES’,  AFIXP  =  ’NO’  means  the  control  function 
is  updated  at  every  adaptive  iteration. 

PARAMETER 

INTCYL:  Defaulted  value  is  0  and  represents  the  number  of  time  steps  at 
which  the  adaptation  is  performed. 

PARAMETER 

NUMCYL:  Defaulted  value  is  999  and  represents  the  interval  of  time  step 
between  adaptations. 

PARAMETER 

MAXINT:  Defaulted  value  is  9999  and  indicates  the  number  of  time  step  at 
which  the  last  adaptation  is  performed. 

PARAMETER 

QUALITY:  Defaulted  value  is  ’NO’. 

QUALITY  =  ’MEASURE’  means  the  grid  is  being  adapted  to  the  quality 
measure  variables  alone. 

QUALITY  —  ’YES’  means  the  grid  is  being  adapted  to  the  combination  of 
quality  measure  variables  and  flow  solution  variables. 
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EXAMPLE  1 

This  example  shows  a  typical  line  of  the  grid  runstream  for  the  case  where 
the  grid  is  being  adapted  statically  to  one  of  the  quality  measure  variables. 

$  ’INITIAL’,  BLEND=2*’ARC\ALL=’YES\CHECK=’N0\ 

CONTYP=’RADIUS\ITMAX=11,CONFAC=0.1, 

ITMAXA=5, 

AWT=’VARVVAR,,,VAR’, 

CW=1.0,1.0,1.0, 

ASKEW*  1M 

AASPE=0., 

AARCL=0., 

APLAC=0., 

AFIXP='NO’, 

ADAPT=’YES’  $ 

Here  the  grid  is  adapted  to  the  skew  angle  variable,  which  is  indicated  by 
ASKEW=1.0  and  AWT=’VAR7VAR7VAR\  The  weight  coefficients  are  taken 
to  be  1.0  in  all  directions  by  CW=1. 0,1. 0,1.0.  The  control  functions  are  updated 
at  each  adaptive  iteration  by  setting  AFIXP=’NO\  and  the  total  number  of 
adaptations  is  ITMAXA=5.  Adaptation  is  enabled  by  ADAPT  =  ’YES’. 
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EXAMPLE  2 

This  example  shows  some  typical  lines  of  the  grid  runstream  for  the  case 
where  the  grid  is  being  adapted  statically  to  one  of  the  flow  solution  variables. 

$  ’INITIAL’,  BLEND=2*’ARC\ALL=’YES\CHECK=’NO\ 

CONTYP=’RADIUS\ITMAX=11,CONFAC=0.1, 

ITMAXA=3, 

AWT=’VAR7GRAD7CURV\ 

CW=0.2,0.5,0.3, 

RHO=1.,0.,0., 

RHOU=0., 

RHOV=0., 

RHOW=0., 

RHOE=0.,0.,0„ 

AFIXP=’YES\ 

ADAPT=’YES’  $ 

C  THE  CONSTRUCTION  OF  THE  GRID 

$  ’VFILE’  ,  FILNAM  =  ’soln.fmt’  ,  FORM  =  ’LIST’  $ 

$  ’END’  $ 

$  ’ERROR’  $ 

$  ’END’  $ 

Here  the  grid  is  adapted  to  the  density  in  the  1  direction,  to  the  gradient  of 
the  density  in  the  2  direction,  and  the  curvature  of  the  density  in  the  3  direction  of 
the  flow  solution  that  is  indicated  by  RHO=l  and  AWT=’VAR’,’GRAD7CURV’. 
The  weight  coefficients  are  CW=0.2,0.5,0.3.  The  control  functions  are  being 
adapted  from  the  original  geometric  forms  in  this  case  since  AFIXP=’YES’, 
and  the  total  number  of  adaptation  is  ITMAXA=3.  Adaptation  is  enabled  by 
ADAPT  =  ’YES’. 
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EXAMPLE  3 

This  example  demonstates  the  use  of  the  namelist  RESTART  to  generate  a 
restart  file  in  the  adaptive  mode. 

$  ’INITIAL’,  BLEND=2*’ARC’,ALL=’YES\CHECK=’N0\ 

CONTYP= ’RADIUS \ITMAX=11,CONFAC=0.1, 
ITMAXA=3, 

AWT=’VAR\’GRAD’,’CURV’, 

CW=0.2,0.5,0.3, 

RHO=1.,0.,0., 

RHOU=0., 

RHOV=0„ 

RHOW=0., 

RHOE=0..0.,0., 

AFIXP=’YES\ 

ADAPT=’YES’, 

RESTART=’YES’, 

VAROUT=’YES’  $ 

Here  VAROUT  =  ’YES’  means  the  adaptive  variable  has  been  saved  in  the 


restart  file. 
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EXAMPLE  4 


This  example  shows  some  typical  lines  of  the  grid  rurstream  coresponding 
to  that  of  Example  3  through  the  use  of  the  RESTART  command  in  the  adaptive 
mode. 


$  ’INITIAL’,  BLEND=2*’ARC’,ALL=’YES\CHECK=’NO\ 

CONTYP=’RADIUS’,ITMAX=11,CONFAC=0.1, 

ITMAXA=3, 

AWT=’VAR\’GRAD7CURV\ 

CW=0.2,0.5,0.3, 

RHO=1.,0.,0., 

RHOU=0„ 

RHOV=0„ 

RHOW=0., 

RHOE=0.,Q.,0., 

AFIXP=’YES’, 

ADAPT=’YES’, 

$  ’BLOCK’  ,  SIZE=35,21,1$ 

$  ’BLOCK’  ,  SIZE=80,41,1$ 

$  ’RESTART’  ,  FILNAM  =  ’rsfile’  ,  VARIN  =  ’YES’  $ 

Here  VARIN  =  ’YES’  indicates  the  adaptive  variable  has  been  read  in  for 
use  in  the  adaptation. 
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EXAMPLE  5 

This  example  shows  some  typical  lines  of  the  grid  runstream  coresponding 
to  the  use  of  the  weighted  average  in  the  adaptation  to  the  flow  variables. 

$  ’INITIAL’,  BLEND=2*’ARC\ALL=’YES\CHECK=’N0\ 

C0NTYP=’RADIUS\ITMAX=1 1,CONFAC=0. 1 , 
ITMAXA=3, 

AWT=’VAR7GRAD7CURV\ 

CW=0.2,0.5,0.3, 

RHO=l.,l.,l., 

RHOU=0., 

RHOV=0., 

RHOW=0., 

RH0E=1.,1.,1., 

PRES=1.,1.,1„ 

AFIXP=’YES’, 

ADAPT=’YES\ 

Here  the  specification  of  RHO  =  1,1,1,  RHOE  =  1,1,1,  and  PRES  =  1,1,1 
means  the  adaptive  variable  is  computed  as  the  combination  of  density,  energy, 
and  presssure  in  each  direction.  Since  AWT  =  ’VAR’,  ’GRAD’,  ’CURV’,  the  grid 
is  being  adapted  to  this  combination  in  the  1  direction,  to  the  gradient  of  this 
combination  in  the  2  direction,  and  to  the  curvature  of  this  combination  in  the 
3  direction. 
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APPENDIX  C 

INCORPORATION  INTO  FLOW  CODES 

The  procedure  for  coupling  of  the  Adaptive  EAGLE  grid  code  into  any  CFD 
flow  code  can  be  described  as  follows: 

1)  First,  define  all  variables  involved  in  the  adaptive  grid  generation  such  as 
AWT,  RESIN,  QUALITY,  AFIXP,  FACTOR,  CW,  RHO,  RHOU,  RHOV, 
RHOW,  RHOE,  PRES,  VOMA,  ITMAX,  ITMAXA,  INTCYL,  NUMCYL, 
KFILE,  ALPHA,  BETA,  ASKEW,  AASPE,  AARCL,  APLAC,  MAXINT, 
ISTAT,  etc.  at  the  declaration  setion  of  the  main  flow  code.  Definition  of 
these  parameters  can  be  found  in  Appendix  B. 

2)  Set  up  a  namelist  to  read  in  all  of  those  variables  defined  in  1. 

3)  Set  up  two  “if”  statements  inside  the  time  step  loop  of  the  main  flow  code, 
one  to  check  at  what  time  step  the  last  adaptation  is  performed  and  the 
other  to  check  how  often  the  adaptation  is  performed. 

Inside  these  “if”  statements 

4)  Call  subroutine  SETVAR. 

5)  Call  subroutine  ELLGEN. 

6)  Call  subroutine  RAGRID. 

The  functions  of  these  subroutines  are  described  below: 

ROUTINE  CALLED  BY  FUNCTION 

SETVAR  Main  flow  code  Computes  the  adaptive  variable  based  on  the  cur¬ 

rent  flow  solution.  After  computing  the  adaptive 
variable  this  subroutine  then  writes  out  these  val¬ 
ues  into  scratch  file  KFILE. 

WRTVAR  SETVAR  Writes  out  the  adaptive  variable  into  scratch  file 

KFILE. 
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RAGRID 

Main  flow  code 

Reads  in  the  adaptive  variable  from  scratch  file 

KFILE. 

ELLGEN 

Main  flow  code 

Generates  a  new  adapted  grid  based  on  flow  vari¬ 
ables  in  the  Dynamic  Adaptation  process. 

There  are  several  subroutines  that  are  called  by  ELLGEN  and  their  functions 

are  described  below: 

ROUTINE 

CALLED  BY 

FUNCTION 

CCDRA & 

CCDWA 

SSD 

When  the  namelist  KSTORE  =  ’CORE’  &  ADAPT 

=  ’YES’,  these  subroutines  are  called  and  have  the 

same  function  as  CCDR  &  CCDW  in  the  original 

case. 

JACBCK 

ELLGEN 

Checks  for  a  twisted  grid. 

QUAL2D  & 

QUAL3D 

ELLGEN 

Computes  the  quality  measures  of  the  grid. 

REDVAR 

ELLGEN 

Reads  in  the  adaptive  variable  that  has  been  writ¬ 
ten  out  on  scratch  file  KFILE  for  use  in  the  adap¬ 
tive  process. 

SSDRA  & 

SSDWA 

SSD 

When  the  namelist  KSTORE  =  ’CORE’  &  ADAPT 

=  ’YES’,  these  subroutines  are  called  and  have  the 

same  function  as  SSDR  &  SSDW  in  the  original 

case. 

REDVAR 

ELLGEN 

Called  when  the  static  adaptation  to  the  available 

flow  variables  is  performed.  It  reads  in  the  flow 

variables  file. 

REDRES 


ELLGEN 


Reads  in  the  restart  file  ’rsfile’. 
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SETIMP 

ELLGEN 

Sets  control  function  values  on  boundaries. 

SETIMR 

ELLGEN 

Sets  coordinates  values  at  Neumann,  image,  and 

reflective  points. 

SETIMV 

ELLGEN 

Sets  the  adaptive  variable  values  at  image  points 

equal  to  the  current  values  at  the  corresponding 

object  points,  sets  weight  values  at  Neumann  and 

reflective  boundaries  and  special  points,  and  also 

extrapolates  to  other  boundary  points. 

SETIMW 

ELLGEN 

Has  the  same  function  as  SETIMV  but  with  the 

weight  values  instead. 

SMOOTHW 

ELLGEN 

After  the  weight  values  are  set  on  the  field  and 

boundaries,  this  subroutine  smooths  these  weight 

values. 

STOREP 

ELLGEN 

When  KSTORE  =  ’FILE’  &  ADAPT  =  ’YES’, 

it  is  called  to  store  the  original  control  function 

values  in  array  for  later  use. 

STOREPW 

ELLGEN 

This  subroutine  has  the  same  function  as  STOREP 

for  KSTORE  =  ’CORE’  k  ADAPT  =  ’YES’  . 

WEIT2D  & 

WEIT3D 

ELLGEN 

Computes  the  weight  function  values  inside  the 

field. 

WEITCON 

ELLGEN 

Performs  the  linear  combination  of  the  original 

control  function  values  that  have  been  saved  by 

STOREP  or  STOREPW  and  the  weight  function 

values. 

WAGRID 

ELLGEN 

Writes  out  the  new  grid  into  scratch  file  KFILE. 
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EXAMPLE 

A  simple  FORTRAN  code  below  demonstrates  the  incorporation  of  the  Adaptive 
EAGLE  grid  into  CFD  flow  codes. 

Begin  of  Adaptive  Section  #1 
(DECLARATION  OF  ALL  PARAMETERS) 

End  of  Adaptive  Section  #1 

DO  10  NTIME=1,NTMAX  (The  time  step  loop  of  main  flow  code) 

Begin  of  Adaptive  Section  #2 

IF(INTCYL  .LE.  MAXINT)  THEN  (Check  for  the  last  adaptation) 
IF(NTIME  .EQ.  INTCYL)  THEN  (Check  for  the  first  adaptation) 
SETVAR  subroutine  computes  the  adaptive  variable  based  on  the  current  flow 
solutions. 

IF( QUALITY  .NE.  ’MEASURE’)  CALL  SETVAR 
ELLGEN  subroutine  performs  the  adaptive  grid  that  bases  on  the  old  grid  and 
the  variables  of  this  INTCYL  time  step. 

CALL  ELLGEN 
RESIN  =  ’NO’ 

RAG  RID  subroutine  reads  in  the  new  grid  that  comes  out  from  ELLGEN. 
CALL  RAGRID 

INTCYL  =  INTCYL  +  NUMCYL 
IF(INTCYL  .GT.  NTMAX)  INTCYL  =  0 
END  IF 
END  IF 

End  of  Adaptive  Section  #2 


10  CONTINUE 
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APPENDIX  D 

FLOW  CODE  ENHANCEMENTS 

An  option  has  been  included  that  allows  the  user  to  select  a  refinement  or 
coarsening  of  the  grid  in  any  of  the  three  coordinate  directions  in  any  block.  The 
grid  selection  is  made  by  the  value  assigned  to  the  character  array  REGRID  in 
the  namelist  FINPUT.  The  following  choices  exist  for  the  array  elements. 

PARAMETER 

REGRID(IC,IB)  =  ’REFINE’.  The  grid  in  block  number  IB  is  to  be  refined 
in  coordinate  direction  IC. 

REGRID(IC,IB)  =  ’COARSE’.  The  grid  in  block  number  IB  is  to  be  coars¬ 
ened  in  coordinate  direction  IC. 

REGRID(IC,IB)  =  ’NO’.  Input  grid  dimensions  are  used  in  coordinate 
direction  IC  of  block  number  IB.  (Default) 

For  example,  REGRID(2,3)  =  ’REFINE’  would  instruct  the  program  to  refine 
block  number  3  in  the  coordinate  direction  2;  that  is,  refine  in  the  J  direction  if 
the  coordinate  directions  are  denoted  I,  J,  and  K. 

The  user  has  the  option  of  creating  a  restart  file  after  the  final  cycle  of  the 
flow  solver  has  been  completed.  When  a  restart  file  is  to  be  written  by  the  flow 
code,  the  solution  may  be  either  output  on  the  actual  computational  grid  as  it 
exists  after  any  grid  refinement  (or  coarsening)  or  interpolated  back  onto  the 
original  grid  that  was  read  into  the  flow  code.  This  choice  is  determined  by  the 
value  of  the  character  variable  DUMP  in  namelist  FINPUT  as  follows. 

PARAMETER 

DUMP  =  ’NEWGRID’.  The  solution  is  output  with  the  grid  that  was  used 
in  the  computations. 
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DUMP  =  ’OLDGRID’.  The  solution  is  output  with  the  grid  that  was  read 
into  the  flow  code.  (Default) 

DUMP  =  ’NO’  means  no  restart  file  is  created. 

Selected  surfaces  or  blocks  of  grid  points  and  solution  values  can  be  written 
out  in  plot3d  formatted  files.  The  character  variable  PL0T3D  in  the  namelist 
FINPUT  is  used  to  indicate  whether  plot3d  data  is  to  be  read  by  the  program. 

PARAMETER 

PL0T3D  =  ’YES’.  PL0T3D  data  is  included  in  the  input  file. 

PLOT3D  =  ’NO’  means  there  is  no  PL0T3D  data  in  the  input  file.  (Default) 

The  block  number  and  the  array  limits  for  each  block  of  plot3d  grid  and  solution 
is  input  into  variables  in  the  namelist  PL0T3D.  The  input  is  terminated  by 
setting  the  block  number  to  zero. 

EXAMPLE 

$  ’PL0T3D’,  BLKA  =  1,  STARTA  =  1,1,1,  ENDA  =  41,41,1  $ 

$  ’PLOT3D’,  BLKA  =  1,  STARTA  =  41,41,1,  ENDA  =  81,41,1  $ 

$  ’PLOT3D’,  BLKA  =  0  $ 

This  data  would  create  a  two-block  plot3d  data  set  from  a  single  surface  grid  in 
the  MISSE  program. 

The  program  can  compute  a  truncation  error  estimate  and  print  out  the 
values  of  the  estimate  by  block  using  various  norms.  The  estimate  can  also  be 
output  in  the  form  of  a  plot3d  solution  file.  The  character  variable  ERROR  in 
the  namelist  FINPUT  has  the  following  options. 

PARAMETER 

ERROR  =  'YES’.  A  truncation  error  estimate  is  computed  and  included  in 
the  output. 

ERROR  =  ’NO’.  A  truncation  error  is  not  computed.  (Default) 


48 


File  names  can  now  be  included  in  the  input  file  for  the  MISSE  code.  This 
option  is  indicated  in  the  value  of  the  character  variable  FNAMES  in  the  namelist 
FINPUT.  The  file  names  are  input  through  the  namelist  FILES. 

PARAMETER 

FNAMES  =  ’YES’  means  the  namelist  FILES  has  values  read  from  the  input 

file. 

FNAMES  =  ’NO’  means  the  default  file  names  are  used.  (Default) 

The  actual  default  values  in  the  MISSE  code  are  used  in  the  following  example 
of  input  for  the  namelist  FILES. 

EXAMPLE 

$  ’FILES’,  GRID  =  ’MISSE.GRD’,  REST  =  ’MISSE.RES’,  P3DGRD  = 
’PL0T3D.GRD’,  P3DSOL  =  ’PLOT3D.SOL’,  P3DERR  =  ’PL0T3D.ERR’  $ 

The  file  GRID  is  the  list  formatted  grid  file  from  the  EAGLE  grid  code.  It 
may  (OUTER  =  ’YES’)  or  may  not  (OUTER  =  ’NO’)  contain  a  surrounding 
layer  of  points  for  each  block.  The  file  REST  is  the  restart  file.  The  remaining 
files  are  plot3d  formatted  files  that  axe  used  to  display  the  solution  and/or  error 
estimates  in  various  regions.  They  are  only  used  if  the  PL0T3D  =  ’YES’  option 
is  exercised. 

A  program  called  COMP  has  been  written  to  compare  numerical  solutions. 
While  it  is  not  a  part  of  the  MISSE  code,  it  is  intended  to  be  used  with  MISSE 
restart  files.  The  program  was  written  to  compare  coarse  and  fine  grid  solutions, 
but  in  fact  can  be  used  to  compare  any  two  solutions  defined  on  or  interpolated 
onto  the  same  grid.  The  input  consists  of  two  MISSE  input  files  and  two  MISSE 
restart  files  for  two  separate  runs  of  the  flow  code.  Along  with  the  input  from 
the  two  flow  codes,  there  are  two  namelist  input  lists,  CINPUT  and  CPLOT, 
which  inputs  information  about  files  and  plotting.  A  sample  namelist  input  file 
is  as  follows: 
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$  ’CINPUT,  P3D  =  ’YES’,  REST1  =  ’MISSE1.RES’,  REST2  =  ’MISSE2.RES’, 
P3DGRD  =  ’COMP.GRD’,  P3DERR  =  ’COMP.ERR’,  DATA1  =  ’MISSE1.DAT’, 
DATA2  =  ’MISSE2.DAT’  $ 

$  ’CPLOT’,  BLKA  =  1,  STARTA  =  1,1,1,  ENDA  =  81,41,1  $ 

$  ’CPLOT’,  BLKA  =  0  $ 

The  input  variables  are  similar  to  the  MISSE  input  and  interpreted  in  the 
obvious  way.  For  example  P3D  =  ’YES’  indicates  that  plot3d  output  is  written 
and  the  other  variables  in  CINPUT  are  the  data  and  restart  files  from  the  MISSE 
solutions.  The  variables  in  the  namelist  CPLOT  are  the  same  as  were  used  in 
the  namelist  PLOT  in  the  MISSE  code.  The  output  consists  of  a  listing  of  the 
parameter  values  used  in  computing  the  two  solutions  and  an  indication  of  when 
there  is  a  difference  in  any  parameter  in  the  two  runs  of  the  MISSE  code.  As  with 
the  flow  solver,  there  is  the  option  for  plot3d  formatted  output  of  the  difference 
in  certain  flow  variables  and  the  complete  solution  vector. 
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Figure  3.2.  Laplacian  with 
V  V  =  0,  V2<f  <  0  . 
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Figure  5.  Contour  plots  of  skew  angle,  aspect  ratio,  Laplacian  of  the  initial  grid. 
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Figvire  6.  Contour  plots  of  the  adaptive  grid,  adaptation  to  skewness 
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Figure  8.  Contour  plots  of  the  adaptive  grid,  adaptation  to  aspect  ratio. 
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Figure  12.  Contour  plots  of  the  adaptive  grid,  adaptation  to  skewness, 
aspect  ratio,  arc  length,  and  Laplacian  of  the  grid. 
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in  the  adaption  to  gradient  and  curvature  of  the  density. 
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Figure  15.  Contour  plots  of  density,  control  functions,  and  weight  functions 
in  the  adaption  to  curvature  and  gradient  of  the  density. 
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Figure  16.  Initial  grid  and  adaptive  grids  for  different  parameters 
and  adaptive  mechanisms. 
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Figure  17.  Contour  plots  of  density  on  initial  and  adaptive  grids. 
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Figure  19.  Pressure  coefficients  of  the  lower  and  upper  surfaces. 
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Figure  20.  Contour  plots  of  density  on  initial  and  adaptive  grids. 


Iptlc  grid  I  fMMt  lv«  grid 


Figure  21.  Contour  plots  of  pressure  on  initial  and  adaptive  grids. 
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Figure  22.  Adaptation  with  AWT=GRAD,GRAD,  PRES=1.1.  CW=0.7,0.7. 
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Convergence  History 
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Figure  23.  Convergence  history  of  the  initial  and  adaptive  grids  solutions. 
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Figure  24.  Contour  plots  of  density  on  fine,  coarse,  and  adaptive  coarse  grids. 
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Figure  25.  Contour  plots  of  the  pressure  on  fine,  coarse,  and  adaptive  grids. 
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Pressure  Coefficient  (Lower  wall) 


Non-adaptation  (Fine  Vs.  Coarse) 


Figure  26.  Pressure  coefficients  of  the  lower  wall  obtained  from 
coarse  and  fine  grids  without  adaptation. 
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Figure  27.  Adaptation  with  AWT=VAR,GRAD.  RHO=U,  CW=0.5.0.5. 
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Figure  28.  Contour  plots  of  density  on  adaptive  grids. 
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Pressure  Coefficient  (Lower  wall) 
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Figure  30.  Adaptation  with  AWT=CURV,CURV,  RHO=l,l,  PRES=1.1  CW=0.7.0. 
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Figure  32.  Adaptation  with  AWT=VAR,GRAD,  RHO=0,1,  PRES=0.1  CW=0  0.0.9. 
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Figure  33.  Convergence  history  of  the  fine,  coarse,  and  adaptive  grid  solutions. 
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Figure  34.  Contour  plots  of  pressure  on  initial  and  adaptive  grids. 
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Figure  35.  Contour  plots  of  error  estimation  on  initial  and  adaptive  grids. 


Figure  39.  Contour  plots  of  pressure  on  a  fine  grid 


Figure  41.  Contour  plots  of  pressure  on  a  coarse  grid, 
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